Europhysics Letters 



PREPRINT 



Dynamics of a Rigid Rod in a Glassy Medium 

Angel J. Moreno 1,2 and Walter Kob 1 

1 Laboratoire des Verres. Universite de Montpellier II. Place E. Bataillon. CC 069. 
F-34095 Montpellier, France. 

2 Dipartimento di Fisica and INFM Udr and Center for Statistical Mechanics and Com- 
plexity. Universita di Roma "La Sapienza". Piazzale Aldo Moro 2. 1-00185 Roma, Italy. 
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Abstract. - We present simulations of the motion of a single rigid rod in a disordered static 
2d-array of disk-like obstacles. The rotational, Dr, and center-of-mass translational, Dcm, 
diffusion constants are calculated for a wide range of rod length L and density of obstacles 
p. It is found that Dcm follows the behavior predicted by kinetic theory for a hard disk with 
an effective radius R(L). A dynamic crossover is observed in Dr. for L comparable to the 
typical distance between neighboring obstacles d nn . Using arguments from kinetic theory and 
reptation, we rationalize the scaling laws, dynamic exponents, and prefactors observed for Dr. 
In analogy with the enhanced translational diffusion observed in deeply supercooled liquids, 
the Stokes-Einstein-Debye relation is violated for L > 0.6d nll . 



Lorentz gases are systems in which a single classical particle moves through a disordered 
array of static obstacles. They are used as simple models for the dynamics of a light particle 
in a disordered environment which presents a much slower dynamics. With increasing density 
of the obstacles, dynamic correlations and memory effects start to become important for the 
motion of the diffusing particle, and the system shows the typical features of the dynamics of 
supercooled liquids, such as a transition to a non-ergodic phase of zero diffusivity [1,2]. The 
theoretical description of the dynamics of the Lorentz gas is highly non-trivial since, e.g., the 
diffusion constant or the correlation functions are non-analytical functions of the density of 
obstacles [1,2]. 

Diffusing particles and obstacles are generally modeled as disks or spheres in two and three 
dimensions, respectively. In this Letter we use molecular dynamics simulations to investigate, 
at low and moderate densities, a generalization of the Lorentz gas, namely a model in which 
the diffusing particle is a rigid rod. In contrast to the usual models, this system has an 
orientational degree of freedom, which for sufficiently long rods is expected to be strongly 
slowed down by steric hindrance. Furthermore this system will allow to gain insight into 
the relaxation dynamics of a nonspherical probe molecule immersed in a simple liquid, an 
experimental technique that is used, e.g., in photobleaching experiments to study the dynamics 
of glass- forming liquids [3], or the investigation of the coupling between the rotational and 
translational degrees of freedom in such liquids [4] . 

In the following we will consider for the array of obstacles a homogeneous "glassy" config- 
uration instead of the more usual random structure. In this way we avoid the effects induced 
by the very different length scales present in a completely random medium [5], and which 
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complicate the physical interpretation of the observed dynamic features. For simplicity, sim- 
ulations have been done in two dimensions, reducing the degrees of freedom to a rotational 
and two translational ones. However, we think that the physical picture proposed here for the 
interpretation of the the obtained results is also valid in a three dimensional system. 

The rigid rod, of mass M, was modeled as TV aligned beads of equal mass m = M/N, with 
a bond length 2a. The rod length is thus given by L — (27V — l)er. In order to obtain the glassy 
host medium, we equilibrated at a density p — 0.77ct~ 2 and at temperature T = e/k^ a two- 
dimensional array of soft disks interacting via a potential V(r) = e(a/r) 12 . This procedure 
produced an homogeneous liquid-like configuration. The latter was then permanently frozen 
and was expanded to obtain the desired density of obstacles, defined as p = n Q b s /l 2 ox , with 
n Q bs the number of obstacles and Ibox the length of the square simulation box used for periodic 
boundary conditions. The same soft-disk potential V(r) was used for the interaction between 
the beads and the obstacles. For computational efficiency, V(r) was truncated and shifted at a 
cutoff distance of 2.5a. In the following, space and time will be measured in the reduced units 
a and (a 2 m/e) 1 ^ 2 , respectively. The rod was equilibrated at T = e/k-Q. After the equilibration, 
a production run was done at constant energy. For statistical average, runs were carried out 
for typically 600-1000 different realizations of the ensemble single rod-obstacles. Runs covered 
10 6 time units, corresponding to (1 — 5) x 10 8 time steps, depending on the step size used for 
the different rod lengths and densities. Run times were significantly longer than the relaxation 
times of the system. 

Fig. 1 shows typical trajectories of the rod center-of-mass (CM) for p = 10~ 2 and differ- 
ent rod lengths. A transition from a typical random walk for short rods to a polygonal-like 
structure for long ones is observed. The latter consists of nearly longitudinal motions inside 
corridors delimited by the neighboring obstacles, connected by vertices corresponding to re- 
gions where the rod changes from a corridor to another one. The very different nature of these 
trajectories suggest a change in the transport mechanism when increasing the rod length at 
constant p. Thus, while for short rods simple rotational and translational Brownian dynamics 
is expected, for long ones rotations and transversal translations are strongly hindered, and 
one has reptation-like motion inside tubes formed by neighboring obstacles and which gives 
rise to the nearly straight long paths observed in the trajectories. 

In order to understand this change in the dynamics we will make use of kinetic theory and 
reptation arguments. Consider first the translational dynamics. We describe the rod by an 
"effective hard rod" formed by a rectangle of width 2a and length L — a and two semi-circles of 
radius a attached at both ends, leading to an effective rod length L + a. This allows us to set 
the size of the obstacles to zero, i.e. to treat them as point particles. If tp is the angle formed 
by the rod axis and the CM velocity, it is easy to see that the projection of the effective rod 
in the transversal direction is 2a + (L — a) sin-0. Integration over all the orientations yields 
an average projection 2R, with R = (L + (n — l)a)/n. The isotropic nature of the Brownian 
regime at low densities suggests to substitute the real system by an equivalent Lorentz gas of 
hard disks having this latter radius. Kinetic theory and mode coupling calculations for hard 
disks of radius R and mean velocity (v) in a 2d- medium with density of point obstacles p lead, 
in first order, to the following result, valid at low and moderate densities [1]: 

D C M/D° CM = l + (32/9w)p*lnp* , (1) 

with the reduced density p* = pR 2 , and D^ M = 3(v)/16pR is the low-density limit, cor- 
responding to Brownian motion -i.e., no correlation between consecutive collisions. Since 
for a Maxwell-Boltzmann distribution the mean CM velocity of the rod is given by (v) = 
(nk^T /2M) 1 / 2 , the corresponding low-density limit for the CM translational diffusion con- 
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stant of the rod, within the picture of equivalent hard disks of radius R = (L + (it — l)a)/ir, 
is: 



and Dcm is obtained from Eqs. (1) and (2) by taking p* = p(L + (it — l)cr) 2 /tt 2 . 

Points in Fig. 2a correspond to the ratio Dqm/Dq M obtained from the simulations, where 
Dqm is calculated as the long-time limit of ((Ar) 2 )/4t, with ((Ar) 2 } the mean squared CM 
displacement. Brackets denote the ensemble average. As can be seen in Fig. 2b, within 
the p— and L— range investigated, the Dcm investigated cover about 4 orders of magnitude. 
Also included in Fig. 2a, is the theoretical prediction given by Eqs. (1) and (2). A good 
agreement with the simulation data is obtained for reduced densities p* smaller than £2 ~ 
0.2, thus showing that in that range, the picture of the equivalent hard disk is correct even 
quantitatively and in particular, the logarithmic corrections to D CM in Eq. (1) give account 
for the observed decrease of -Dcm / D CM at intermediate values of p* . Since the typical distance 
between neighboring obstacles, d nn , is given by d nn = p~ x / 2 , p* = £2 corresponds to a rod 
length L w 1.4<f nn , and to a diameter for the equivalent hard disk of 2R « 0.9d nn . For 
hard disks of this size, one has very strong correlation effects and the first order expansion in 
Eq. (1) is not valid any more [1]. Indeed the continuation of Eq. (1) to larger values of p* 
leads to an unphysical increase of the ratio Dqm/ D CM (see the theoretical curve in Fig. 2a). 
In contrast to this, the increase of this ratio observed in the simulation results does have a 
physical origin. (Note that this increase is found only in the ratio and that Dcm decreases 
with increasing p* , see Fig. 2b.) However, it is not possible to rationalize this increase within 
the picture of equivalent hard disks of radius R= (L + (n — 1)<j)/tt. Instead, reptation of the 
effective rod allows the system to stay ergodic, which would not be the case for hard disks at 
large p*, and hence to increase the ratio beyond 1.0. We will give below evidence for reptation 
as the transport mechanism for large p*. 

For the case of the rotational dynamics, a calculation from first principles, that would 
require to take into account the interactions of the obstacles with all the beads forming the 
rod, is not straightforward. However, it is well-known that heuristic approximations in the 
framework of elementary kinetic theory [6], as the one we are going to discuss now, lead, for 
very low densities, i.e., in the Brownian regime, to correct results within a factor 0(1). For 
a given ensemble of n rods with CM velocity v and angular velocity u), we calculate n(9), 
the number of rods that rotate a given angle 9 without colliding. A differential equation for 
this quantity can be obtained by taking into account that the number of rods undergoing a 
collision between 9 and 9 + d9 is n(9) multiplied by the probability of a collision in d6. By 
taking advantage of the equivalent hard disk picture, such a probability can be obtained as p 
times the area swept by the equivalent hard disk, 2Rdx = 2(L + (it — l)a)dx/n, where dx and 
d9 are related through dx = vd6/u. Thus we obtain the differential equation n(9+d9)—n(9) = 
dn = —2n(9)pRvd9/ui 1 and after integration, n(9) — n cx-p(—2pRv9/uj). The mean angular 
free path between collisions is obtained as (9) = (1/no) / dvf(v) J f(w)dw f™°9dn, where f(v) 
and f(u) are the Maxwell-Boltzmann distributions at temperature T of the CM and angular 
velocities respectively. Integration leads to (9) = (u))(l/v)/(2pR). For a Maxwell-Boltzmann 
distribution we have (u>) — (2k^T '/nl) 1 / 2 , with / the moment of inertia of the rod, and 
(1/v) = (7rM/2fcBT) 1 / 2 . Therefore, from all these results, and the expression from kinetic 
theory for the rotational diffusion constant, Or = (9}(u), we obtain: 
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Dr is determined from the simulation as the long time limit of ((A</>) 2 )/2t, with ((A(/>) 2 ) 
the mean squared angular displacement. In Fig. 3 we show its dependence on p* (left set of 
data points and left ordinate). The thick solid line for p* < £i « 0.04 represents Eq. (3) and 
shows a good agreement with the simulation results. This value of p* = £i corresponds to 
a rod length L « 0.6d nn . As a further test of the validity of this heuristic approach, we can 
calculate the low-density limit D^ M following the same arguments for the equivalent system 
of hard disks, leading to a decay rate for a displacement r, n(r) — no cxp(— 2pRr), a mean 
free path I = (2pR)~ 1 , and a diffusion constant Dqm = l(v)/2= {v)/ApR, i.e., a factor which 
is only 4/3 larger than the rigorous result D^ M = 3(v)/16pR. In Fig. 3 a crossover to a 
different dynamic regime is observed around p* = £i. As can be seen in Fig. 2, this value of 
p* corresponds to a decay of a 15% from the low-density limit of Dqm, again consistent with 
the breakdown of the picture of Brownian dynamics. 

For the regime p* > £i we now follow reptation arguments similar to those developed for 
solutions of rigid linear polymers [7]. Thus, we assume that at a given time, the hindrance 
effects induced by the neighboring obstacles result in the confinement of the rod in the center 
of a tube of width S and length L+a equal to that of the effective rod. We decompose the long- 
time rotational diffusion as a sum of elementary processes i, each one having a persistence time 
t% inside a tube, where the rod sweeps a typical angle The rotational diffusion constant 
is estimated as Dr = fi 2 /2r, where the average values of ft and r are taken. Due to the 
homogeneous nature of the host medium, only a small dispersion is expected for tti and t,, 
and the latter equation should be a reasonable approximation. The persistence time r can 
be estimated as the time that the rod takes to cover its effective length L + a. Once the rod 
has moved this distance, the old tube has been left and a new tube is defined by the new 
neighboring obstacles. Therefore, r can be obtained as r = (L + a) 2 / 4D^ M , with D^ M the 
low-density limit given by Eq. (2), since we are referring to the translational dynamics of the 
rod inside the free space of the tube. On the other hand, the angle swept inside the tube can 
be estimated as ft — 25/ (L + a). Due to the homogeneous configuration of the host medium, 
the width S can be estimated as S — d nn — 2a = p~ x l 2 — 2a, where the term 2a is an excluded 
area effect due to the finite width of the effective rod. From all these results, Dr. = tt 2 /2r 
follows the scaling law: 



A good agreement with the simulation data is obtained, as shown by the thick solid line 
in Fig. 3 for p* > £i w 0.04. Only a small, though systematic, deviation is observed at 
the largest values of p* , suggesting that in this range the model is somewhat oversimplified. 
Except for these largest values, the agreement is particularly good for p* > £2 ~ 0.2, i.e., 
for the range where reptation was suggested as mechanism of ergodicity restoring for the 
translational motion. For £1 < p* < £2, corresponding to a rod length 0.6c? nn < L < 1.4d nn , 
reptation also gives a reasonable description of the simulation data, although a more detailed 
theory might be required to describe this crossover regime. 

It must be stressed that there is no functional continuity between Eqs. (3) and (4) and 
hence a more rigorous model should be introduced to obtain a unified description at both 
sides of the crossover point p* = £1. Nevertheless, we find that the curves for D R predicted 
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by (3) and (4) actually cross at a value p* = £ co in agreement with the value of £i « 0.04 
obtained from the simulations. It is straightforward to see that this common value is given 
by £ co = tt- 2 {1 + (tt - 1)(<j/L)} 2 {2((t/L) + (2M/3nI) 1 / 2 (l + (a / L)) 2 L}~ 2 . This expression 
can be simplified by taking into account that, except for the highest investigated densities, 
the data in the reptation regime fulfill the condition a -C L and therefore one can make the 
approximation I = ML 2 / '12. Using these two approximations we obtain for £ co the result 
p* = (87r) _1 w 0.04, in agreement with the observed value of £i. From all these results, at this 
point it can be concluded that Brownian motion and reptation are the transport mechanisms 
respectively below and above p* = £i, the crossover taking place when the rod length is 
comparable to the typical distance between neighbouring obstacles. It is also noteworthy of 
remark that a consistent description has been achieved without any fit parameter. 

Some insight about the coupling between the rotational and translational degrees of free- 
dom can be obtained by investigating the validity of the Stokes-Einstein-Debye (SED) rela- 
tion [8] DcmTr = c, where c is a constant and tr is the rotational relaxation time. We define 
this latter quantity as the time at which the correlation function (cosA0(i)) has decayed to 
1/e. In the limit of small p* , where collisions are infrequent, one expects that the rod performs 
several full rotations between two consecutive collisions. In this case, the physical mechanism 
for the angular decorrelation is not the collisions but the free rotation of the rod between them. 
Therefore, in this limit tr can be calculated from the relation cos((cj)tr) = 1/e, leading to 
tr = (7T//2/CB7 1 ) 1 / 2 cos _1 (l/e). A good agreement of this latter equation with the simulation 
data has been obtained for p* < £i, confirming the validity of the hypothesis (see Ref. [5]). 
Then, together with the low-density limit of Dcm given by Eq. (2), we obtain the following 
scaling law for the SED-relation: 

D° cm t rP [L + (tt - l)a] (Mj If I 2 = (3^ 2 /32) cos-^l/e) . (5) 

As shown in Fig. 4, the simulation data follow nicely this equation for p* < £i, i.e., for L < 
0.6d nn . Data for £i < p* < £ 2 are systematically below the constant value (37r 2 /32) cos _1 (l/e), 
showing that the SED-relation is no more fulfilled in this range of p* (the reason for this is the 
presence of the logarithmic correction given in Eq. (1)). Finally, a crossover to a power law 
~ (p*) 2 ' 4 is obtained at p* = £2- (Note that presently we do not have an explanation for this 
power-law.) Thus, the SED-relation is clearly violated for rods that are significantly longer 
than d nn , in analogy with the phenomenon of enhanced translational diffusion observed in 
supercooled liquids -with temperature as control parameter-, close to the glass transition [9], 
where the value of c increases with decreasing temperature, and that is explained by the differ- 
ent averaging of the rotational and translational degrees of freedom [10] . This interpretation is 
clearly evidenced by the trajectory of the center-of-mass for L = 59 in Fig. 1 and the picture 
of reptation at large values of p* , where due to steric hindrance, the rod is able to perform 
only small rotations while it moves along nearly straight long paths. In this way, reptation 
leads to an enhancement of the translational decorrelation in comparison with the rotational 
one. 

In summary, we have shown that arguments from kinetic theory and reptation are able 
to give a very good acount for the translational and rotational dynamics of the considered 
generalization of the Lorentz gas. In particular the analysis show that in this system the 
violation of the Stokes-Einstein-Debye relation is due to the breakdown of the isotropy of the 
rotational dynamics, a result that will be useful for the interpretation of the experimental 
findings of this violation. 
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Fig. 1 - Typical trajectories of the rod CM at p — 10 -2 , for L — 3, 39 and 59. The simulation time is 
isim = 10 6 , with a time interval of 500 between consecutive plotted points. The arrow indicates the 
length scale. Note that for the sake of clarity the obstacles are not shown. 
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Fig. 2 - (a): Ratio between the CM translational diffusion constant Dcm obtained from the sim- 
ulations and the low-density limit Dq M for the equivalent hard disk (symbols). The solid curve 
corresponds to the theoretical prediction of kinetic theory for the equivalent hard disk at low and 
moderate densities, see Eq. (1). (b): CM translational diffusion constant vs. the reduced density. 
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Fig. 3 - Scaling laws for Dr. A double-?/ representation is shown for clarity. Points are simulation 
data. Thick solid lines correspond to the theoretical predictions for Brownian and reptation dynamics 
(see text). 
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Fig. 4 - Scaling law for the product Dcmtr and violation of the SED relation at large p* . The dashed 
line marks the theoretical constant value (37r 2 /32) cos _1 (l/e) (see text). The solid line in the large-p* 
regime is a fit to a power law. 
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